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ABSTRACT 

Gamma-ray loud X-ray binaries are binary systems that show non-thermal broadband emission from radio to 
gamma rays. If the system comprises a massive star and a young non-accreting pulsar, their winds will collide 
producing broadband non-thermal emission, most likely originated in the shocked pulsar wind. Thermal X-ray 
emission is expected from the shocked stellar wind, but until now it has neither been detected nor studied in 
the context of gamma-ray binaries. We present a semi-analytic model of the thermal X-ray emission from 
the shocked stellar wind in pulsar gamma-ray binaries, and find that the thermal X-ray emission increases 
monotonically with the pulsar spin-down luminosity, reaching luminosities of the order of 10 33 erg s . The lack 
of thermal features in the X-ray spectrum of gamma-ray binaries can then be used to constrain the properties of 
the pulsar and stellar winds. By fitting the observed X-ray spectra of gamma-ray binaries with a source model 
composed of an absorbed non-thermal power law and the computed thermal X-ray emission, we are able to 
derive upper limits on the spin-down luminosity of the putative pulsar. We applied this method to LS 5039, 
the only gamma-ray binary with a radial, powerful wind, and obtain an upper limit on the pulsar spin-down 
luminosity of ~ 6 x 10 36 erg s _1 . Given the energetic constraints from its high-energy gamma-ray emission, a 
non-thermal to spin-down luminosity ratio very close to unity may be required. 

Subject headings: gamma rays: stars — stars: individual (LS 5039) — stars: winds, mass-loss — X-rays: 
binaries 



1. INTRODUCTION 

In the past few years, gamma-ray binaries have emerged 
as a new category of very high-energy (VHE, E > 100 GeV) 
gamma-ray sources. In most of the cases, they comprise a stel- 
lar mass compact object orbiting a young, massive star. They 
generally differ from typical X-ray binaries in that they exhibit 
luminosities at gamma-ray energies at levels similar or above 
their X-ray luminosity. Five such sources have been already 
detected: LS5039 (Aharonian et al. 2005b), PSRB 1259-63 
(Aharonian et al. 2005a), LSI +61 303 (Albert et al. 2006), 
HESS J0632+057 (Skilton et al. 2009; Bongiorno et al. 201 1; 
Moldon et al. 2011), and 1FGL J1018.6-5856 (Corbet et al. 
201 1; Pavlov et al. 201 1). Evidence of VHE flaring emission 
has been found in Cygnus X-l (Albert et al. 2007), but at 
present lacks confirmation. 

Whereas the pulsar nature of the compact object in 
PSRB 1259-63 is determined thanks to the detection of ra- 
dio pulses, this is not possible in the case of LS 5039 and 
LS I +61 303, in which the radio pulses would be likely free- 
free absorbed in the stellar wind. Broadband radiation model- 
ing in the microquasar scenario has been able to reproduce the 
spectral energy distribution of the source (e.g., Paredes et al. 
2006; Bosch-Ramon et al. 2006), but the lack of accretion 
signatures, and the similarities with PSRB 1259-63, prompted 
other authors to consider the pulsar wind shock (PWS) sce- 
nario for these sources (e.g., Martocchia et al. 2005; Dubus 
2006; Chernyakova et al. 2006). 

The emission from the non-thermal particle population ac- 
celerated in the PWS scenario has been widely studied in 
the radio, X-ray, and high energy (HE) and VHE gamma-ray 
bands (Maraschi & Treves 1981; Tavani & Arons 1997; Kirk 
et al. 1999; Dubus 2006; Khangulyan et al. 2007; Neronov & 
Chernyakova 2007; Bosch-Ramon 2011). Furthermore, the 
possibility of detecting the free pulsar wind from its spec- 



tral signature has also been discussed (e.g., Kirk et al. 1999; 
Khangulyan et al. 2007; Cerutti et al. 2009). However, one of 
the sources of emission that could provide valuable information 
on the characteristics of gamma-ray binaries has been hitherto 
overlooked: the thermal X-ray emission from the shocked stel- 
lar wind. Massive stellar binaries (of types O+O or O+WR) 
are known to produce powerful thermal X-ray emission in the 
wind collision region (WCR), where the supersonic winds of 
the two young stars interact (see, e.g., Pittard et al. 2005, for a 
review). Even though it has been previously considered that 
the detection of X-ray emission lines would be an indication of 
accretion in gamma-ray binaries, similar lines can arise from 
the thermal X-ray emission of the shocked stellar wind in the 
WCR. The complexity of the stellar winds from the Be com- 
panion stars in LSI +61 303 and PSRB 1259-63 makes the 
analysis of the thermal emission from the shocked stellar wind 
difficult. The presence of both a polar wind and a decretion 
disk adds a number of unknown free parameters. On the other 
hand, the radial wind from the O star in LS 5039, also more 
powerful than the one from Be stars, turns this source into 
an ideal candidate to study the observational impact of the 
thermal wind emission. 

In this paper we present a model for the thermal X-rays 
from pulsar gamma-ray binaries inspired on semi-analytical 
models of the X-ray emission from massive binary systems 
(e.g., Antokhin et al. 2004; Parkin & Pittard 2008). This model, 
which allows the study of the dependence of the thermal X- 
ray emission for different stellar and pulsar wind properties, 
will be applied to the binary LS 5039 in order to constrain the 
spin-down luminosity of its putative pulsar. In Section 2, we 
characterize the shocked stellar wind region, approximated 
here as the contact discontinuity between the pulsar and the 
stellar winds, and present a model of its thermal X-ray emis- 
sion in Section 3. In Section 4 we apply the model to LS 5039 
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and present spin-down luminosity upper limits in Section 5. 
Finally, we discuss the implications of the results in Section 6. 

2. THE DYNAMICAL MODEL 

In the PWS scenario, the pulsar wind and the stellar wind 
collide forming an interaction region bounded on either side 
by the reverse shocks of the winds. The shocked winds are 
separated by a surface called contact discontinuity (CD), lo- 
cated where the ram pressure components perpendicular to 
this surface are in equilibrium: = p p± . Both winds are 
assumed here to be radial. 

The ram pressure of the relativistic pulsar wind can be 
represented at any point at a distance r p from the pulsar as 
Pp( r p) = L s a/(4ncrp), where L s( j is the spin-down luminosity 
of the pulsar. In the case of the stellar wind, we will assume a 
/3-velocity law in the radial direction, v(r*) = Voo(l - R+/r+)P, 
where is the distance to the center of the star, 7?* is the 
radius of the star and is the terminal velocity of the stellar 
wind. The azimuthal component of the wind will be neglected. 
We have considered /3 = 0.8, which has been shown to be a 
good approximation of the outer structure of the winds of hot 
massive stars (Pauldrach et al. 1986). Assuming a constant 
mass-loss rate from the star M, the density of the stellar wind 
is p — M 7(4717-^(7*)). The ram pressure of the stellar wind is 
given by p*(r*) = Mv(r*)/(47irJ). 

The balance of perpendicular components of the ram pres- 
sure that defines the position of the CD 
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can be rearranged to obtain the dimensionless parameter rj 
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where and 8 p are the angles between the CD surface and 
the directions toward the center of the star and the pulsar, re- 
spectively. The coordinates x and y correspond to the position 
along the line of centers and the distance from it, respectively 
(see Figure 1). Following Antokhin et al. (2004), we can obtain 
the differential equation that describes the shape of the CD as 
a function of r] from Equation (2), 



d.v 
dV 



Drl(x,y) yjj](x,y) 
rl(x,y) ^r](x,y) + r 2 p (x,y) 



(3) 



where x(y) is the function describing the shape of the CD and 
D is the binary separation distance. The boundary condition 
(dx/dy | v ^o — > 0) is xq = D/(l + ^jrj), the location of pressure 
equilibrium along the symmetry axis, also known as stagnation 
radius. For the case of constant wind velocities, rj will be con- 
stant, and xq and x(y) can be computed analytically. However, 
given the y6-velocity law we adopted for the stellar wind, a 
numerical integration is required to obtain the shape of the CD 
since the stellar wind may not yet have reached its terminal 
velocity at the contact surface. Figure 1 shows the resulting 
shape of the CD in comparison to the star-pulsar system as a 
function of t]^ = (Z, sd /c)/(Mv M ). 

At the scales of interest, and given the very high speed of 
the stellar wind of massive stars, the role of orbital motion 
in the shape of the CD is not expected to be significant. The 
aberration or skew angle p of the shock cap is defined as the 
angle between the axis of symmetry of the shock cap and the 




Figure 1. Illustration of the dynamical model for a close pulsar gamma-ray 
binary. The ~10R o star (but not the pulsar) is at scale. The thick black line 
indicates the shape of the CD for rja, = 0.04. The line CP is tangent to the CD 
at the point P, for which the angles and radii used in the derivation of the CD 
are shown. The gray lines show the shape of the CD for values of ^„ spaced 
logarithmically in the range 0.003-0.3. 



line of centers of the star and compact object. This angle is 
given by tan^u = v or b/v w , where v or b is the orbital velocity and 
v w is the speed of the slowest wind (i.e., the wind of the star) 
at the contact surface. For the orbital parameters of LS 5039, 
for example, the skew angle is below 23° at all moments along 
the orbit. 

3. THERMAL X-RAY EMISSION 

3.1. Cooling regime of the shocked stellar wind 

On either side of the CD, the stellar and pulsar winds develop 
shocks. The shocked pulsar wind is the candidate location for 
particle acceleration, which would give rise to the observed 
broadband non-thermal emission (Tavani & Arons 1997). On 
the other side of the CD, the stellar wind develops a strong 
shock where the wind material is slowed down and heated. 
The cooling efficiency of the gas determines the width of the 
cooling layer. On one extreme, if radiative cooling is very 
efficient, the gas will be rapidly cooled after the shock and will 
collapse into a thin dense layer on the CD. The rapid cooling 
assures that all the wind kinetic energy flux perpendicular to 
the shock is transformed into thermal emission. On the other 
extreme, if radiative cooling is slow enough, the hot gas will 
flow along the CD and cool adiabatically through expansion. 
In the adiabatic regime, a smaller fraction of the kinetic energy 
of the stellar wind is emitted in the X-ray band. Stevens et al. 
(1992) proposed the parameter^ = f ra dAesc as a measure of 
the cooling regime of the shocked material in colliding wind 
binaries, x values significantly below (above) unity indicate 
that the shocked stellar wind cools radiatively (adiabatically), 
whereas \ values around unity indicate that a fraction of the 
incoming kinetic energy will be radiated away and the rest will 
be dissipated through adiabatical expansion. 

3.2. Estimation of the expected X-Ray luminosity 

Before we enter into the details of the calculation of the ther- 
mal X-ray spectrum and luminosity, we will roughly derive the 
expected X-ray luminosity as a function of the parameters of 
the system. As mentioned above, we expect the X-ray luminos- 
ity to be a fraction / of the wind kinetic luminosity crossing 
the shock: Lx = fL^ n - The value of / is strongly dependent 
on the value of x at the shock, and as a first approximation can 
be taken as / » 1/(1 + x) ( see Equation (7)). Additionally, 
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the kinetic luminosity deposited on the shock, L^?, will be a 
fraction of the total kinetic luminosity of the stellar wind, L??, 
determined by the fraction of the total solid angle of the stellar 
wind subtended by the shock. For a cone with a half-opening 
angle of 8, we obtain Lx ~ -cos ff)L^, which can be sim- 
plified up to OiO 4 ) as L x « f^6 2 L^. To obtain the estimate 
of the luminosity, we will here consider that the wind reaches 
the shock front at terminal velocity. The angle 6 can then be 
derived from Equation (3) as 7tt]/(1 + rj), which on a first-order 
approximation for small values of rj (i.e., dominance of the 
stellar wind) will be 9 « nrj. Using these approximations we 
obtain Lx ~ f\n 2 rfL}^ n . As seen in Equation (2), 77 is propor- 
tional to L s( j. Therefore, we expect the emitted thermal X-ray 
luminosity from the shocked stellar wind to behave roughly 
as: 



1.2 x 10 



32 



^sd 



10 36 ergs~ 



M 



10- 7 M G yr- 



erg s 



(4) 



where we have taken the wind terminal velocity as v M 
2400 km s . and assumed that half of L. sh is converted into 

kin 

X-ray luminosity, i.e., / = 0.5. This latter assumption is consis- 
tent with a cooling regime between the adiabatic and radiative 
regimes, i.e.,^- w 1. 

The above approximations will only hold for small values of 
77, i.e., for the case of a dominant stellar wind. When the ram 
pressure of the pulsar wind is similar to that of the stellar wind 
and 77 approaches unity, the half-opening angle of the CD will 
approach n/2. At this limit, the X-ray luminosity will be even 
higher that for the previous case, and of the order of 



3 x 10 34 



M 



10- 7 M yr- 



erg s 



(5) 



Gamma-ray binaries are non-thermal X-ray sources with 
luminosities of the order of several 10 33 ergs _1 . Therefore, 
as Equations (4) and (5) show, pulsar spin-down luminosities 
moderately higher than 10 36 erg s" 1 would imply a significant 
thermal component in their X-ray spectrum. 

3.3. X-Ray emission model details 

For a strong, steady, standing shock with a compression ratio 
Po/Pw = Vwj_/voj_ = 4 (where the subindices w and indicate 
the properties of the gas immediately before and after the 
shock, respectively), the gas will have an immediate postshock 
temperature of 



kT = — im p v. 



,2 

W-L 



1.21 



0.62 



Vwi. 



1000 km s" 



keV, (6) 



where v WJ _ indicates the preshock speed of the wind component 
perpendicular to the shock front and yU is the average atomic 
weight. A value of jj — 0.62 is appropriate for a fully ionized 
medium with solar abundances (Anders & Grevesse 1989). 

The X-ray emission of the shocked stellar gas is computed 
assuming that, in the emitting plasma, conditions are equal 
to those at the immediate postshock position, which is true 
as long as the gas remains subsonic. This will provide an 
estimate of how much of the kinetic energy of the incoming 
wind is radiated away. Assuming the shocked gas acts as a 
calorimeter, in which all the incoming perpendicular kinetic 
energy is either radiated away as X-ray emission or dissipated 
through adiabatic expansion, the fraction of the kinetic energy 



converted into X-ray emission will be 



f rad + f < 



r sh 

-1 Man' 



(7) 



where f la d is the radiative cooling timescale, f esc is the escape 
or flow dynamics timescale and L^? is the incoming kinetic 
energy. Given that E cc t , the fraction of times may be 
interpreted as the fraction of the radiated energy over the total 
energy lost owing to radiation and escape. 

For a gas cooling with a volumetric emission rate A(T), the 
radiative cooling timescale is given approximately by 



frad - 



kT n 



n A(T )' 



(8) 



where no is the postshock number density, and in a strong shock 
will be four times the preshock wind number density. Instead 
of using a pre-calculated table for A(T), it can be calculated 
using the mekal code used to compute the emitted spectrum 
(see below). 

The escape timescale is the time the heated gas takes to flow 
away from the postshock location. Assuming that the stellar 
wind overpowers the pulsar wind, we can use the distance from 
the CD to the pulsar as the radius of curvature of the shock 
front, and thus as a typical flow dynamics distance scale. The 
escape timescale can then be approximated as f esc = r p /v . 
We assume that the wind velocity component parallel to the 
shock is unaffected after traversing it, so that the immediate 
postshock speed may be found as = (v w± /4) 2 + vjL. 

We used the mekal code as implemented in the spex software 
package (Kaastra et al. 1996) to generate a look-up table of 
emission spectra for optically thin hot gas in collisional ioniza- 
tion equilibrium as a function of the gas temperature 70 We 
populated the look-up table with emission spectra in the 0.05- 
15 keV range, taking 250 logarithmically spaced temperatures 
from 10 4 to 10 9 K. Additionally, for each of the temperatures, 
the volumetric emissivity rate A(T) was computed by inte- 
grating the emission between photon energies of 0.001 and 
50keV. 

In order to compute the total cumulative emission from the 
shocked stellar wind we assumed that the X-ray emitter is 
infinitely thin and located on the CD defined through Equa- 
tion (3). We defined a grid of 2000 bins in the y-direction, away 
from the line of centers, and 30 azimuthal bins on the shock 
front. For each segment of the grid the immediate postshock 
temperature To and X-ray luminosity was obtained through 
Eqs. (6) and (7). Using the X-ray emission look-up table, 
we obtained the intrinsic X-ray emission spectrum from each 
segment. Finally, we computed the neutral hydrogen column 
density owing to the non-shocked stellar wind, from each of 
the positions of the segments on the shock front to the ob- 
server. We assumed that neither the enhanced density of the 
shocked stellar wind nor the free pulsar wind contribute sig- 
nificantly to the total column density. The column density of 
the former can be estimated to be two orders of magnitude 
below the wind column density (Szostek & Dubus 201 1), and 
the latter is too rarefied to absorb X-rays. We then computed 
the energy-dependent attenuation factor corresponding to the 
found column density through the wabs subroutine provided 
with Xspec vl2. In addition, the emission from those areal 
segments eclipsed by the star was discarded. The absorbed 
emission from each of the areal segments was added together to 
obtain the cumulative emitted spectrum of the shocked stellar 
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wind. 

From the calculations described in this section, we can see 
what the general characteristics of the shocked stellar wind 
emission will be. Along the line of centers of the system we 
will find both the highest temperatures and the highest densi- 
ties owing to the perpendicular impact of the wind on the shock 
front. The hardest X-rays will come from this location, but its 
volume is relatively small and its luminosity may not dominate 
the final emission. Away from this point temperatures will drop 
rapidly as the incoming wind will not be perpendicular to the 
shock, whereas the volume increases, and density decreases. 
This region will produce a softer X-ray spectrum. The final 
spectrum will be a so-called multi-temperature X-ray spec- 
trum, resulting from the emission of the stratified temperature 
structure of the shocked gas. 

3.4. Caveats of the model 

The semi-analytic model detailed in the previous sections is 
robust given that the assumptions made are fairly conservative. 
However, there are some factors that affect the WCR that can- 
not be accounted for in our semi-analytic framework. Detailed 
hydrodynamical simulations would be required to take into 
account effects such as the following. 

• In the radiative cooling regime the CD may be disrupted 
by thin-shell instabilities (Stevens et al. 1992; Pittard 
2009). Additionally, stellar wind clumping would likely 
disturb the shocked wind region to some extent. If this 
perturbation is propagated to the shocked pulsar wind, it 
would have implications for both the thermal and non- 
thermal emission. 

• We modeled the radiative output of the postshock gas 
through the ratio of radiative to total cooling timescales 
in Equation (7), and assumed an infinitely thin cooling 
layer. Both the intrinsic X-ray luminosity and cooling 
layer width could be calculated self-consistently through 
hydrodynamical simulations. Whereas the results are 
likely to be very similar in the radiative cooling regime 
limit, they will probably differ in the adiabatic limit. In 
the adiabatic limit of the postshock flow, the flow may 
become supersonic at a distance smaller than r p , cooling 
the gas and reducing the radiative outcome. This effect 
would not be important at scales close to the apex, but 
could affect the contribution produced at the periphery 
of the system. 

3.5. Method of spin-down luminosity upper limit derivation 

To accurately compare the output of the model described in 
Sections 2 and 3 with the observed X-ray spectra of gamma- 
ray binaries, the computed X-ray spectra can be formatted as a 
FITS table model in order to load it into the spectral analysis 
software package Sherpa (Freeman et al. 2001). The FITS 
table model is generated from the cumulative spectra for 25 
pulsar spin-down luminosities spaced logarithmically between 
10 36 erg s~'and 5 x 10 37 erg s _1 , taking a constant value of the 
mass-loss rate and orbital inclination. 

The X-ray spectra observed at a certain orbital phase can 
be assigned a source model corresponding to an interstellar 
medium neutral hydrogen absorption model (A), affecting both 
the non-thermal power law (P) and the computed shocked 
stellar wind thermal emission (Th): Aism X(P + Th). The goal 
is then to find for which value of the spin-down luminosity the 
thermal component begins to affect the shape of the spectrum 



and the fit worsens significantly. We used the confidence 
Sherpa command to evaluate the 99.7% (3<x) confidence level 
range for acceptable values of the spin-down luminosity. If 
no thermal features are present in the observed spectrum, the 
lower bound will be consistent with zero, i.e., lack of thermal 
component. The upper bound, on the other hand, provides us 
with the upper limit to the spin-down luminosity of the pulsar. 

4. APPLICATION OF THE MODEL TO LS 5039 

The system LS 5039 is located at 2.5 + 0.1 kpc and contains a 
compact object with a mass between 1.4 and 5 M Q , orbiting an 
06.5V((f)) donor star every 3.90603 ± 0.00017 days (Casares 
et al. 2005) in a mildly eccentric orbit (<? = 0.24 + 0.08; Sarty 
et al. 2011). The detection of elongated asymmetric emission 
in high-resolution radio images was interpreted as mildly rel- 
ativistic ejections from a microquasar jet and prompted its 
identification with an EGRET HE gamma-ray source (Pare- 
des et al. 2000, 2002). However, recent Very Long Baseline 
Array observations by Ribo et al. (2008) show morphological 
changes on short timescales that might be consistent with a 
pulsar binary scenario. It shows variable periodic emission 
in the X-ray, and HE and VHE gamma-ray bands. The X-ray 
lightcurve has an orbital modulation (Takahashi et al. 2009), 
superposed by short timescale features that are quite stable over 
the years (Kishishita et al. 2009). The maximum and minimum 
emission phases in the periodic X-ray lightcurve are located 
around the inferior and superior conjunctions, so a geometrical 
effect for the modulation has been suggested. However, infe- 
rior and superior conjunctions are located close to apastron 
and periastron, and the broadness of the peaks, well fitted by 
a sinusoidal function, makes it difficult to discern between a 
physical or geometric origin of the modulation. Whereas the 
VHE gamma-ray lightcurve has a modulation similar to the 
X-ray one (Aharonian et al. 2006), the HE gamma-ray flux is 
lower at the inferior conjunction (0.45 < (f> < 0.9) and higher at 
superior conjunction (0.9 < <p < 0.45). In addition, at superior 
conjunction it shows a spectral cutoff at E = 1.9 + 0.5 GeV, 
not statistically required for inferior conjunction data (Abdo 
et al. 2009). 

The radial and powerful wind of LS 5039 makes it an ideal 
system to apply the model presented in this work. The shocked 
gas cooling regime (see Section 3.1) will be closer to the 
radiative than to the adiabatic limit, with^ values around unity. 
Therefore, it seems safe to apply the thermal emission model 
presented in Section 3.3 to LS 5039. 

4. 1 . Stellar mass-loss rate and pulsar spin-down luminosity 

in LS 5039 

As explained in Section 2, the shape of the CD depends only 
on the parameter 77 M , which in turn depends only on the ram 
pressure of the winds on either side of the CD. Figure 2 shows 
the shape of the CD in comparison to the orbit of LS 5039 at 
periastron and apastron for different values of rj^. Casares 
et al. (2005) derived that the 06.5V((f)) star in LS 5039 has a 
stellar radius of /?* = 9.3 Rq, its wind has a terminal velocity 
of Voa = 2400 km s~' and the orbital inclination is in the range 
20°-70°. Therefore, the only parameters still unknown that 
affect the shape of the CD, and thus the thermal X-ray emission, 
are the stellar mass-loss rate, M, and the pulsar spin-down 
luminosity, L sc j. 

Sarty et al. (2011) performed an optical photometric and 
spectroscopic campaign on LS 5039, and determined the mass- 
loss rate of the star to be 3.7 to 4.8 x 10" 7 M yr~' from Ha 
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Figure 2. Sketch of the orbit of the compact object around the companion star 
in LS 5039. The star (but not the pulsar) is at scale. The two marked positions 
correspond to periastron (P) and apastron (A). For each of these positions the 
shape of the CD is shown for tj m = 0.0025, 0.025 and 0.25 (blue, red, and 
green lines, respectively). 



line measurements. However, Ho- line fitting is a diagnostic 
known to be affected by wind clumping, and owing to its p 2 
dependence could easily overpredict mass-loss rates by a factor 
~2 (Markova et al. 2004). Bosch-Ramon (2010) studied the 
effect of wind clumping in the absorption of X-ray emission 
from LS 5039, but it must be noted that Sarty et al. (201 1) did 
not detect any evidence of dense wind clumps in the analysis 
of stellar wind emission lines. The lack of variability of the 
column density measured in X-ray observations along the orbit 
may also be used to constrain the maximum value of mass-loss 
rate. Bosch-Ramon et al. (2007) concluded that, considering 
a non-thermal X-ray emitter inside the binary system, the 
mass-loss rate could not exceed a few 10~ 8 M yr _1 . Adopting 
a model in which the non-thermal emission is extended, as 
expected in the PWS scenario, this upper limit is relaxed up to 
1.5 X 10" 7 Mo yr _1 (Szostek & Dubus 2011). In the following 
calculations we considered the whole range of mass-loss rates 
discussed here, i.e., from 5xlO -8 M Q yr~' to 4.8xl0~ 7 M yr _1 . 
For some cases, in which a single representative mass-loss rate 
is used to test other properties of the thermal emission, we 
chose the midpoint of the range: M = 2.65 x 10~ 7 M Q yr _1 . 

The spin-down luminosity of the putative pulsar in LS 5039 
is not known as it has never been measured directly. A lower 
limit may be derived from efficiency considerations taking into 
account the HE gamma-ray emission detected by Fermi/LAT 
(Abdo et al. 2009). If the HE gamma-ray emission has a 
magnetospheric origin, the detected flux above 100 Me V of 
2.55 x 10~ 10 erg cirT 2 s _1 implies a spin-down luminosity be- 
tween 0.2 and 30 x 10 36 erg s _1 , depending on the beaming 
model, with higher values favored by the high energy spectral 
cutoff (see Zabalza et al. 201 1, for an explanation of the deriva- 
tion of these values for the case of LS 1+61 303). However, the 
observed orbital modulation of the HE gamma-ray emission, 
and the lack of pulsations, indicate that its origin might not be 
magnetospheric but rather IC. 



The minimum non-thermal luminosity required for the Fermi 
emission of LS 5039 can also be estimated in the context of IC 
radiation. To obtain a conservative estimate, one can assume 
that electrons are injected only in the relevant energy range, ~1- 
100 GeV, their velocities are isotropically distributed, and have 
time to cool down only through IC. Stellar IC leading to pho- 
tons in the Fermi band occurs in the Thomson regime. Since 
this process is strongly anisotropic (Lie °c (1 — cos #tc) (p+1 ^ 2 
with p being the electron power-law index; e.g., Dermer et al. 
1992), one needs to account for the IC interaction angle (f?ic) 
between the observer line of sight and the seed photon di- 
rection. Under these conditions, the injection luminosity of 
electrons can be written as 



2{a + 1)(1 -cos6 lc ) a ' 

where the factor 1.7 accounts for the luminosity radiated below 
100 MeV, a = (p+l)/2, and 6ic < n/2 forphases around the in- 
ferior conjunction. Taking Lc e v ~ 2x 10 35 erg s _1 during those 
orbital phases (Abdo et al. 2009), p ~ 3 to explain the Fermi 
spectrum, and a rather conservative 8\c ~ 60°, one obtains 
L e ~ 2 x 10 36 ergs~'. However, escape, adiabatic, and per- 
haps synchrotron losses could easily increase L e by a factor of 
several. In addition, not all the L s d will go to non-thermal par- 
ticles 1 , and the injection electron distribution is expected to be 
broader than -1-100 GeV. From all this, L sd > 2 x 10 37 erg s -1 
seems more realistic. Doppler boosting has not been accounted 
for but, although it can reduce the required energy budget in 
certain orbital phases, it will increase it in others. Note that for 
6>ic « 90° and phases between inferior and superior conjunc- 
tion, i.e., with the flow bulk motion unlikely pointing to us, 
still L e > 10 36 ergs _1 , and thus L sd > 10 37 ergs^. The GeV 
emitter may be extended, thus relaxing the f?ic dependence, but 
it should not change significantly our conclusion. 

4.2. Comparison to X-Ray data 

LS 5039 has been observed at different positions along its 
orbit by many X-ray observatories (see Zabalza et al. 2008, for 
a summary). Recently, a very long observation with Suzaku 
provided an uninterrupted lightcurve of the source during an 
orbit and a half (Takahashi et al. 2009), and a comparison with 
previous observations indicated that it was quite stable over 
periods of up to a decade (Kishishita et al. 2009). The Suzaku 
observation provides a very good orbital coverage, but we have 
chosen, to compare the output of our numerical model, two 
XMM-Newton observations performed in 2005 (Bosch-Ramon 
et al. 2007). The superior effective area of XMM-Newton gives 
us the possibility to extract very good soft X-ray spectra in rel- 
atively short observations of duration A(p ~ 0.04. These short 
observations assure that the ambient conditions do not vary 
significantly and the measured spectra can be accurately com- 
pared to spectra calculated for a given phase. These two ~15 ks 
observations were performed at periastron (Obsld 0202950301, 
orbital phases 0.02-0.05) and apastron (Obsld 0202950201, 
orbital phases 0.49-0.53) during the same orbital period, thus 
probing the two extreme separations between the star and the 
compact object. 

For each of the two observations we processed the data 
through the standard procedure with version 10.0 of the XMM- 
Newton Science Analysis Software. We filtered out the periods 

1 Not only may the acceleration efficiency be smaller than 100%, but for 
higher spin-down luminosities the opening angle of the CD will increase and 
only around half of the pulsar wind luminosity will be shocked. 
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Table 1 

Results of the X-Ray Spectral Analysis 



Orb. Phase 


Parameter 


Source Model 






Aism X P 


A KM x(P + Th) 


Periastron 


r 

F{Pf 
F(Thf 


1.57 ±0.04 
7.1 ±0.6 


1.56 ±0.03 
6.6 ± 0.4 
0.37 


Apastron 


r 

F(Pf 
F{Thf 


I. 47 ±0.03 

II. 4 ±0.7 


1.42 ±0.03 
10.9 ±0.6 
0.58 




N H b 
L s d 

x 2 /doi 


0.67 ± 0.02 
467.5/612 


0.69 ± 0.02 
4.09 x 10 36 ergs-' 
477.1/612 



11 Unabsorbed X-ray fluxes of the non-thermal power law (P) 
and the thermal shocked stellar wind (TK) components in the 
range 0.3-10keV are given in units of 10~ 12 erg cirT 2 s -1 . 
b Interstellar column density is in units of 10 22 cirT 2 . 



of high particle background for each of the three detectors 
aboard XMM-Newton, and extracted source and background 
spectra from the clean event files, as well as the corresponding 
redistribution matrix and ancillary response files. We binned 
the extracted source spectra to a minimum of 20 counts per bin, 
avoiding an oversampling of the detector energy resolution by 
more than a factor three. 

Before introducing our model of the thermal X-ray emission 
of the shocked stellar wind, we performed a fit of the data to 
an absorbed power law using the spectral analysis package 
Sherpa bundled with CIAO 4.3. The fit results for this strictly 
non-thermal source model may be found in the third column 
of Table 1. 

We assigned each of the two groups of spectra (periastron (p) 
and apastron (a)) the source model described in Sections 3.5: 
Aism x (P p ' d + Th p ' d ). Whereas the normalization and photon 
index parameters of the power-law component were left free 
for the spectral fitting, the interstellar absorption and the spin- 
down luminosity of the pulsar were kept constant between the 
periastron and apastron source models. We note that such a 
choice for the source model implies the assumption that the 
non-thermal emission is not significantly affected by circum- 
stellar absorption, whereas it affects the thermal emission and 
is already computed in the component Th. This would happen 
for either extended non-thermal emitters or low circumstellar 
absorption (see Section 6. 1 below for a discussion). 

As expected from the previous estimate (see Section 3.2), the 
X-ray luminosity of the thermal component increased mono- 
tonically with the value of the spin-down luminosity of the 
pulsar. In order to derive the spin-down luminosity upper limit, 
we performed a simultaneous fit of the six spectra (three de- 
tectors for each of the two orbital phases). In all cases, the 
best fit corresponded to the lowest level of emission from the 
thermal component, i.e., the non-thermal component alone pro- 
vided the best fit to the data. The upper bound, as explained in 
Sec. 3.5, corresponds to the upper limit to the pulsar spin-down 
luminosity. 

5. RESULTS 

5.1. Thermal X-Ray spectra 

We used the dynamical, radiative, and absorption models 
described in Sections. 2 and 3 to calculate the thermal X-ray 
spectra emitted by the shocked stellar wind of LS 5039. For 
the properties of this system we find that the cooling regime 



„ io 34 

au 



Apastron 


\ K k K \ \ ^ \ \ \ \ \ * 

\ \\ \\ X \X \\ ; 


- — Periastron 




Estimate 






^^^^^ 




\\ Excluded by\\^ 










■ i ■ , , , 





10 36 1 37 

Pulsar spin-down luminosity L s( j [ergs -1 ] 



Figure 3. Thermal emission luminosity in the 0.3-10 keV range as a function 
of the spin-down luminosity of the pulsar. The fluxes at periastron (solid) and 
apastron (dashed) are shown. In addition, the rough estimate of Equation (4) 
is shown as a gray line. The range of pulsar spin-down luminosities excluded 
by the thermal emission is shown as a hatched region. In all cases, a mass-loss 
rate of 2.65 x 10~ 7 M Q yr 1 was assumed. 

of the shocked gas is an intermediate regime, with x around 
unity. As an example, for M = 2.65 x 10~ 7 M Q yr ~' and 
L s d = 3x 10 36 erg s _1 (^ M = 0.025), we obtain radiative cooling 
and escape timescales at the CD apex during periastron of 9.1 x 
10 3 s and 7.7 x 10 3 s, respectively, resulting my — 1.18. The 
maximum postshock temperature reached is Tq - 2.43 keV 
at the apex, with a postshock number density of «o - 2 x 
10 10 cnT 3 . 

In accordance with the estimates of Section 3.2, we found 
that in general the computed thermal X-ray luminosities in- 
crease monotonically with the spin-down luminosity of the 
pulsar. In Figure 3 we show the 0.3-10 keV X-ray luminos- 
ity of the thermal emission from the shocked stellar wind 
as a function of the pulsar spin-down luminosity, taking 
M = 2.65 x 10 _7 M o yr _1 . For scenarios with t]^ <k 1 the 
behavior is generally in good agreement with Equation (4), 
and the thermal X-ray luminosity can be acceptably explained 
by a function of the form L x °c L" d with a - 1.381 ± 0.004. 
This value is significantly lower than that predicted by the 
estimate of Equation (4). This could be related to the lack of 
correction for the obliquity of the incoming flow with respect 
to the shock in the estimate, which could result in an under- 
estimation for lower 7700 and an overestimation at higher 1}^, 
where a larger fraction of the incoming wind impacts obliquely. 
At higher rj^, the proximity of the CD to the stellar surface 
at periastron provokes a decrease in the velocity of the wind 
crossing the shock, and therefore a decrease in the total X-ray 
luminosity. 

Figure 4 shows the obtained intrinsic and absorbed spectra 
for different values of the spin-down luminosity from 3 X 
10 35 ergs~' to 3 x 10 37 ergs _1 . The mass-loss rate of 2.65 x 
10~ 7 M Q yr assumed for these spectra provides a dense matter 
field in the binary. Its effect can be seen as a very strong 
absorption at energies below 1 keV for the periastron spectrum. 
The column density at this phase is of the order of 10 22 cm -2 , 
even higher that the measured interstellar one. In contrast, the 
apastron spectrum for the highest spin-down luminosity shows 
negligible circumstellar absorption. At apastron, the pulsar 
would be close to inferior conjunction, and the high opening 
angle of the CD would mean that most of the thermal X-ray 
emission only has to travel through the unshocked pulsar wind 
to reach the observer (see Figure 2), thus being barely absorbed. 
As discussed in Section 3.3, the hardest X-rays will come from 
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Figure 4. X-ray emission computed from the shocked stellar wind of LS 5039 
at periastron (top) and apastron (bottom), taking M = 2.65 X 10 -7 M Q yr~' . 
The different spectra for each orbital phase correspond to different values of 
pulsar spin-down luminosity L sl j: 0.3 (blue, bottom), 3 (red, middle), and 30 
(green, top) times 10 36 ergs~', corresponding to = 0.0025, 0.025, and 
0.25, respectively. The thin dashed lines indicate the intrinsic spectra, whereas 
the thick solid ones are the ones absorbed owing to the stellar wind. 

the region around the apex of the CD. For high values of the 
spin-down luminosity, the apex of the CD at periastron will be 
close to the stellar surface, and the velocity of the incoming 
wind will be low owing to the wind /3-velocity law. Therefore, 
the X-ray emission from the apex will be softer, an effect that 
can be seen in the high energy part of the topmost spectrum of 
Figure 4. For the parameters used to compute this spectrum, 
the stagnation point of the CD is at only 0.39 stellar radii from 
the surface of the star, and the wind velocity at the shock is 
about one third of the terminal velocity. Note how this effect is 
not apparent in the corresponding apastron spectra, for which 
the stagnation point is much farther from the stellar surface. 

5.2. Spin-down luminosity upper limits 

To derive the upper limits of the pulsar spin-down luminosity, 
we chose to explore the stellar mass-loss rate values mentioned 
in Section 4.1, corresponding to the estimates of Bosch-Ramon 
et al. (2007), Szostek & Dubus (2011), and Sarty et al. (2011): 
5 x 10~ 8 Mq yr 1 , 1 .5 x 10~ 7 M yr 1 , and 4.25 x 10~ 7 M yr 1 , 
respectively. The orbital inclination (z = 20°-70°; Casares et al. 
2005) will affect the importance of circumstellar absorption 
at periastron in the emitted spectrum. The result of the pulsar 
spin-down luminosity upper limit determination, for the whole 
ranges of mass-loss rate and orbital inclination, is shown in 
Figure 5. A summary of the results for the extremes and 
midpoints of the ranges can be seen in Table 2. 

We found that pulsar spin-down luminosities above 6 x 
10 36 ergs 1 are excluded by the lack of spectral thermal fea- 




20° 30° 40° 50° 60° 70° 

Orbital inclination 



Figure 5. Spin-down luminosity upper limits as a function of orbital inclina- 
tion and stellar mass-loss rate. 



Table 2 

Pulsar Spin-down Luminosity Upper Limits 



Mass-loss Rate 


Inclination 


L s( j Upper Limit 


(M Q yr-') 




(10 36 ergs-') 




20° 


3.3 


5 x 10~ 8 


45° 


4.3 




70° 


5.6 




20° 


4.5 


1.5 x 10~ 7 


45° 


4.9 




70° 


4.8 




20° 


4.4 


4.25 x 10~ 7 


45° 


4.0 




70° 


3.9 



tures in the observed periastron and apastron spectra, for any 
combination of mass-loss rate and orbital inclination. This 
upper limit corresponds to the extreme inclination of 70°, and 
for moderate orbital inclinations the upper limit is at all times 
below 5.5 x 10 36 erg s _1 . Although the wind kinetic luminosity 
available for M — 5 x 10~ 8 M Q yr is lower than for the other 
cases, the upper limits derived are at least as strict as for higher 
mass-loss rates. We found that the reduced absorption owing to 
the lower density of the stellar wind results in brighter emission 
below 1 ke V. In the last column of Table 1 , we show the best-fit 
results of the thermal plus non-thermal source model when fix- 
ing the spin-down luminosity at the upper limit we found. The 
thermal component was computed taking the midpoint values 
of mass-loss rate and orbital inclination. As compared to the 
non-thermal source model, the best-fit power law is slightly 
harder and its unabsorbed flux is ~5% lower. An illustration 
of the thermal and non-thermal contributions to the observed 
spectrum is shown in Figure 6. 

6. DISCUSSION 

6.1. On the stellar mass-loss rate 

The stellar mass-loss rate plays a double role in the PWS 
scenario. On one hand it determines the shape of the CD and 
the kinetic luminosity that is shocked onto the CD. On the other 
hand, it determines the importance of photoelectric absorption 
in the X-ray spectra of the thermal and non-thermal radiation 
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Energy [keV] 

Figure 6. X-ray spectra of LS 5039 at periastron (top) and apastron (bottom). 
The open circles are the data from the pn detector. The MOS data are fit 
simultaneously but not shown here for clarity. The lines indicate the best- 
fit spectral model for the spin-down luminosity fixed at the upper limit for 
moderate values of orbital inclination (i = 45°) and stellar mass-loss rate 
(M = 2.65 X 10~ 7 M o yr'): L sd = 4.09 X 10 36 ergs"'. The thick solid line 
is the sum of the non-thermal power-law component (thin dashed) and the 
shocked stellar wind thermal component (thin solid). The energy resolution of 
the thermal component is lower than in Figure 4, which results from folding 
the spectral models with the energy redistribution matrix of the detector. 



sources. 

The lack of X-ray absorption variability in LS 5039 along 
the orbit indicates that, if the non-thermal emitter is inside the 
orbital system, the stellar wind density must be low. However, 
a very low value will result in the pulsar wind overpowering 
the stellar wind for reasonable values of the pulsar spin-down 
luminosity. The maximum mass-loss rate compatible with a 
variation of the neutral hydrogen column density along the 
orbit below 10 21 cirT 2 (Bosch-Ramon et al. 2007) is of ~ 5 x 
10~ 8 M yr _1 . For such a low mass-loss rate, the pulsar wind 
could impact on the stellar surface at periastron for reasonable 
values of spin-down luminosity. However, the lack of orbital 
variations of the UV spectra of LS 5039 (McSwain et al. 2004) 
indicate that the stellar surface is not perturbed at any point 
along the orbit. This places an upper limit on the spin-down 
luminosity of ~ 7 x 10 36 ergs _1 . However, if the star has a 
strong enough surface magnetic field, then magnetic pressure 
could be enough to balance the pulsar wind away from the 
surface (Harding & Gaisser 1990). For a surface magnetic 
field of B+ — 50 G and assuming a dipole field, the spin-down 
luminosity could be as high as 3 x 10 37 ergs~' without the 
CD collapsing onto the stellar surface. Note, however, that 
M = 5 x 10~ 8 M yr _1 and fi* = 50 G imply a wind magnetic 
to kinetic energy ratio close to unity. 

If a significant fraction of the non-thermal emitter is not 



located deep inside the binary system, the photoelectric ab- 
sorption of the thermal and non-thermal X-ray emission could 
then be decoupled (as we have done in Section 3.5), since 
their locations would be different. For this scenario, the stellar 
mass-loss rate could be as high as a few times 10~ 7 M Q yr , as 
measured by Sarty et al. (201 1). Additionally, such a location 
would favor the acceleration of particles up to the very high 
energies required to account for the VHE gamma-ray spectrum, 
as well as avoid strong pair creation absorption of VHE gamma 
rays (Khangulyan et al. 2008; Bosch-Ramon et al. 2008). 

6.2. On the pulsar spin-down luminosity 

The found upper limits of the pulsar spin-down luminosity 
pose strong constraints on further modeling of the PWS sce- 
nario for non-thermal emission. As discussed in Section 4.1, 
the pulsar spin-down luminosity should be at least a few times 
10 36 ergs _1 in order to account for the GeV emission. The 
upper limit we have found from the thermal X-ray emission 
(L s d < 6 x 10 36 ergs _1 ) indicates that the spin-down lumi- 
nosity must be very close to this value in order for the PWS 
scenario to be consistent. A detailed modeling of the GeV 
emission from LS 5039 is required to accurately check the 
consistency of the observed GeV luminosity and the PWS sce- 
nario. Unfortunately, the structure and location of the GeV 
emitter is not well known, and IC emission is very sensitive 
to the 6ic -value, which makes difficult to obtain conclusive 
results. Nevertheless, the lower limit derived here for L e , in 
the IC scenario, is still robust enough to show that, likely, 
L e /L s( j — > 1. It might appear that such a strong constraint can 
only be achieved at very low Tj m , but the effect of the orbital 
motion on the large-scale structure allows for the pulsar wind 
to be terminated even in the direction opposite to the star up 
to rjcn ~ 1 (Bosch-Ramon & Barkov 201 1). We note that the 
constraint on the efficiency is slightly less restrictive for high 
orbital inclinations (i > 45°) and low stellar mass-loss rates 
(M<2x lO^Moyr 1 ). 

The recent detection of PSRB 1259-63 at GeV energies 
(Tarn et al. 201 1; Abdo et al. 201 1) provides a perfect compar- 
ison of the efficiency of spin-down luminosity to non-thermal 
emission conversion efficiency. During its highest state of 
emission the total gamma-ray luminosity was comparable to 
the pulsar spin-down luminosity, which combined with the re- 
cent reclassification of the optical star (Negueruela et al. 201 1) 
results in nearly unbearable energy requirements (Khangulyan 
et al. 201 1). Even considering Doppler boosting, the efficiency 
of the shocked pulsar wind in converting the spin-down lumi- 
nosity into non-thermal HE gamma-ray emission must be quite 
high, and would only weaken the constraint for certain phases 
(see Section 4.1). 

The acceleration and emission processes responsible for the 
high efficiency in the conversion from spin-down luminos- 
ity to HE gamma-ray emission are currently unknown. IC 
emission from particles accelerated beyond the light cylinder 
in the striped pulsar wind is one of the proposed emission 
processes (Petri & Dubus 2011). However, when applied to 
PSRB 1259-63, this process can only account for the low- 
flux state before periastron, and is unable to explain the high 
luminosity detected after the periastron passage without the 
invocation of addition seed photon sources apart from the stel- 
lar companion. The extremely high efficiency required for 
LS 5039, similar to the high-flux state of PSR B 1259-63, in- 
dicates that this process might not be responsible for its HE 
gamma-ray emission. 
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6.3. Concluding remarks 

The semi-analytic model presented here has allowed us to 
gain valuable insight into the nature of LS 5039. The assump- 
tions described in Section 3 are chosen such that the resulting 
level of thermal X-ray emission is fairly conservative and 
thus the upper limit on the spin-down luminosity robust. We 
have shown that the study of the thermal emission from the 
shocked stellar wind can be a powerful diagnostic tool for 
pulsar gamma-ray binaries. Even for cases where no thermal 
features are observed in the X-ray spectrum, it can be used to 
place constraints on important properties of the system such 
as stellar mass-loss rate and pulsar spin-down luminosity. The 
application of such a model to LS 5039 allowed us to conclude 
that, if hosting a non-accreting pulsar, this system has a very ef- 
ficient non-thermal mechanism, with L s( j ~ (3-6) x 10 36 erg s _1 . 



We acknowledge support by the Spanish Ministerio de 
Ciencia e Innovation (MICINN) under grants AYA2010- 
21782-C03-01 and FPA2010-22056-C06-02. V.Z. was sup- 
ported by the Spanish Ministerio de Education through FPU 
grant AP2006-00077. V.Z. and V.B.-R. thank the Max 
Planck Institut fur Kernphysik for its kind hospitality and 
support. The research leading to these results has received 
funding from the European Union Seventh Framework Pro- 
gramme (FP7/2007-2013) under grant agreement PIEF-GA- 
2009-252463. J.M.P. acknowledges financial support from 
ICREA Academia. 

Facilities: XMM (EPIC) 



REFERENCES 



Abdo, A. A., Ackermann, M., Ajello, M., et al. 2009, ApJ, 706, L56 
Abdo, A. A., Ackermann, M., Ajello, M., et al. 201 1, ApJ, 736, LI 1 
Aharonian, F., Akhperjanian, A. G., Aye, K., et al. 2005a, A&A, 442, 1 
Aharonian, F, Akhperjanian, A. G., Aye, K., et al. 2005b, Science, 309, 746 
Aharonian, F, Akhperjanian, A. G., Bazer-Bachi, A. R., et al. 2006, A&A, 
460, 743 

Albert, J., Aliu, E., Anderhub, H., et al. 2006, Science, 312, 1771 
Albert, J., Aliu, E., Anderhub, H., et al. 2007, ApJ, 665, L51 
Anders, E., & Grevesse, N. 1989, Geochim. Cosmochim. Acta, 53, 197 
Antokhin, 1. 1., Owocki, S. P., & Brown, J. C. 2004, ApJ, 611, 434 
Bongiorno, S. D., Falcone, A. D., Stroh, M., et al. 201 1, ApJ, 737, LI 1 



Bosch-Ramon, V. 2010, in ASP Corjf. Ser., Vol. 422, High Energy 
Phenomena in Massive Stars, ed. J. Marti, P. L. Luque-Escamilla, & 
J. A. Combi (San Francisco, CA: ASP), 77 
Bosch-Ramon, V. 201 1, Proc. 25th Texas Symposium (arXiv:l 103.2996) 
Bosch-Ramon, V., & Barkov, M. V. 201 1, A&A, submitted (arXiv:l 105.6236) 
Bosch-Ramon, V., Khangulyan, D., & Aharonian, F. A. 2008, A&A, 489, L21 
Bosch-Ramon, V., Motch, C., Ribo, M., et al. 2007, A&A, 473, 545 
Bosch-Ramon, V., Paredes, J. M., Romero, G. E., & Ribo, M. 2006, A&A, 
459, L25 

Casares, J., Ribo, M., Ribas, I., et al. 2005, MNRAS, 364, 899 
Cerutti, B., Dubus, G., & Henri, G. 2009, A&A, 507, 1217 
Chernyakova, M., Neronov, A., & Walter, R. 2006, MNRAS, 372, 1585 
Corbet, R. H. D., Cheung, C. C, Kerr, M., et al. 201 1, ATel, 3221 
Dermer, C. D., Schlickeiser, R., & Mastichiadis, A. 1992, A&A, 256, L27 
Dubus, G. 2006, A&A, 456, 801 

Freeman, P., Doe, S., & Siemiginowska, A. 2001, Proc. SPIE, 4477, 76 

Harding, A. K., & Gaisser, T. K. 1990, ApJ, 358, 561 

Kaastra, J. S., Mewe, R., & Nieuwenhuijzen, H. 1996, in UV and X-ray 

Spectroscopy of Astrophysical and Laboratory Plasmas, ed. K. Yamashita 

&T. Watanabe, 411^14 
Khangulyan, D., Aharonian, F, Bogovalov, S., & Ribo, M. 201 1, ApJ, 

submitted (arXiv:1104.0211) 
Khangulyan, D., Aharonian, F, & Bosch-Ramon, V. 2008, MNRAS, 383, 467 
Khangulyan, D., Hnatic, S., Aharonian, F, & Bogovalov, S. 2007, MNRAS, 

380, 320 

Kirk, J. G., Ball, L., & Skjaeraasen, O. 1999, Astropart. Phys., 10, 31 
Kishishita, T., Tanaka, T., Uchiyama, Y., & Takahashi, T. 2009, ApJ, 697, LI 
Maraschi, L., & Treves, A. 1981, MNRAS, 194, IP 
Markova, N., Puis, J., Repolust, T., & Markov, H. 2004, A&A, 413, 693 
Martocchia, A., Motch, C, & Negueruela, I. 2005, A&A, 430, 245 
McSwain, M. V., Gies, D. R., Huang, W., et al. 2004, ApJ, 600, 927 
Moldon, J., Ribo, M., & Paredes, J. M. 201 1, ATel, 3180 

Negueruela, I., Ribo, M., Herrero, A., et al. 201 1, ApJ, 732, LI 1 
Neronov, A., & Chernyakova, M. 2007, Ap&SS, 309, 253 
Paredes, J. M., Bosch-Ramon, V., & Romero, G. E. 2006, A&A, 451, 259 
Paredes, J. M., Marti, J., Ribo, M., & Massi, M. 2000, Science, 288, 2340 
Paredes, J. M., Ribo, M., Ros, E., Marti, J., & Massi, M. 2002, A&A, 393, 
L99 

Parkin, E. R., & Pittard, J. M. 2008, MNRAS, 388, 1047 
Pauldrach, A., Puis, J., & Kudritzki, R. P. 1986, A&A, 164, 86 
Pavlov, G. G., Misanovic, Z., Kargaltsev, O., & Garmire, G. P. 201 1, ATel, 
3228 

Petri, J., & Dubus, G. 2011, MNRAS, 1193 
Pittard, J. M. 2009, MNRAS, 396, 1743 

Pittard, J. M., Dougherty, S. M., Coker, R. R, & Corcoran, M. F. 2005, in 
X-Ray and Radio Connections, ed. L. O. Sjouwerman & K. K. Dyer 

Ribo, M., Paredes, J. M., Moldon, J., Marti, J., & Massi, M. 2008, A&A, 481, 
17 

Sarty, G. E., Szalai, T, Kiss, L. L., et al. 2011, MNRAS, 411, 1293 
Skilton, J. L., Pandey-Pommier, M., Hinton, J. A., et al. 2009, MNRAS, 399, 
317 

Stevens, I. R., Blondin, J. M., & Pollock, A. M. T. 1992, ApJ, 386, 265 
Szostek, A., & Dubus, G. 2011, MNRAS, 411, 193 
Takahashi, T, Kishishita, T, Uchiyama, Y, et al. 2009, ApJ, 697, 592 
Tam, P. H. T, Huang, R. H. H, Takata, J., et al. 201 1, ApJ, 736, L10 
Tavani, M., & Arons, J. 1997, ApJ, 477, 439 

Zabalza, V., Paredes, J. M., & Bosch-Ramon, V. 2008, Int. J. Mod. Phys. D, 
17, 1867 

Zabalza, V., Paredes, J. M., & Bosch-Ramon, V. 201 1, A&A, 527, A9 



